// ================================================================================================
// Intersection test for a triangle surface
// ================================================================================================
inline float TriangleSurface::SurfaceIntersection(const Ray& ray, float tmax) {
    Vector4 b = triangle.v2 - triangle.v0;
    Vector4 c = triangle.v1 - triangle.v0;
    Vector4 normal = CrossProduct(c,b);
    float h[3];

	float tplane = -DotProduct((ray.origin-triangle.v0),normal) / DotProduct(ray.direction,normal);
	if (tplane < 0.0f || tplane > tmax) 
        return 0.0f;

    float absx = fabs(normal.x);
    float absy = fabs(normal.y);
    float absz = fabs(normal.z);
    uint k;
    if (absx > absy)
        if (absx > absz) 
            k = 0; /* X */ 
        else 
            k=2; /* Z */
    else
        if (absy > absz) 
            k = 1; /* Y */ 
        else k=2; /* Z */

    uint u = (k+1) % 3; 
    uint v = (k+2) % 3;

	// calc hitpoint
	h[u] = ray.origin[u] + tplane * ray.direction[u];
	h[v] = ray.origin[v] + tplane * ray.direction[v];
    float hu = h[u] - triangle.v0[u];
    float hv = h[v] - triangle.v0[v];

    float denom = (b[u] * c[v]) - (b[v] * c[u]);
	float beta = ((b[u] * hv) - (b[v] * hu)) / denom;
	if (beta < 0.0f) 
        return 0.0f;

	float gamma = ((c[v] * hu) - (c[u] * hv)) / denom;
	if (gamma < 0.0f) 
        return 0.0f;

	if (beta+gamma > 1.0f) 
        return 0.0f;

    return tplane;
}

